#include "../kissfft/kiss_fft.h"
#include "../kissfft/kiss_fftr.h"
//https://github.com/mborgerding/kissfft.git

void generate_triangular_wave(float *samples, int num_samples,  float sample_rate) {
    // const int num_samples = 4096;
    const float frequency = 1000.0f; // 1kHz
    
    // 相位增量 = 频率 / 采样率
    float phase_increment = frequency / sample_rate;
    float phase = 0.0f;

    for (int i = 0; i < num_samples; i++) {
        // 计算当前相位对应的幅值
        float value;
        if (phase < 0.5f) {
            // 上升阶段：0 → 1
            value = 4.0f * phase - 1.0f;
        } else {
            // 下降阶段：1 → -1
            value = 3.0f - 4.0f * phase;
        }
        
        samples[i] = value;
        
        // 更新相位（保证在[0,1)范围内）
        phase += phase_increment;
        if (phase >= 1.0f) {
            phase -= 1.0f;
        }
    }
}
void generate_sine_wave(float *samples, int num_samples, float sample_rate) {
    // const int num_samples = 4096;
    const float frequency = 1000.0f; // 1kHz
    const float two_pi = 2.0f * (float)M_PI;

    for (int i = 0; i < num_samples; i++) {
        // 直接计算每个采样点的相位，避免累计误差
        float t = (float)i / sample_rate;
        samples[i] = sinf(two_pi * frequency * t);
    }
}

int main2()
{
    int N = 1024;
    kiss_fft_cfg cfg = kiss_fft_alloc(N, 0, NULL, NULL);
    kiss_fft_cpx *in = (kiss_fft_cpx*)malloc(N * sizeof(kiss_fft_cpx));
    kiss_fft_cpx *out = (kiss_fft_cpx*)malloc(N * sizeof(kiss_fft_cpx));
    
    float sample[1024];    
    // generate_triangular_wave(sample, 1024, 16000);
    generate_sine_wave(sample, 1024, 16000);

    for (int i = 0; i < N; i++) {
        // in[i].r = (kiss_fft_scalar)((rand() + 10) % 256) / 256.0;
        in[i].r = sample[i];
        in[i].i = 0;
    }
    // 填充输入数据到in数组（需转换为复数格式）
    kiss_fft(cfg, in, out); // 
    
    
    for (int i = 0; i < N; i++) {
        if (out[i].r >= 0.001 || out[i].i >= 0.001) {

        printf("%d  (%4.4f+%4.4fi), \n",i, out[i].r, out[i].i);
        }
    }

    // 处理输出（计算幅值）
    free(in);
    free(out);
    kiss_fft_free(cfg);
    return 0;
}

#include <stdio.h>
static void write_file(const char *name, void *data, int len)
{
    FILE *bin_fp = fopen(name, "wb");
    if (!bin_fp) {
        perror("无法创建二进制文件");
        exit(1);
    }
    // 写入频谱数据（float数组）
    fwrite(data, len, 1, bin_fp);

    fclose(bin_fp);

}
#include <stdint.h>
#include <stdlib.h>
static void write_float(const char *name, float *p, int len)
{
    int16_t *s = (int16_t *)malloc(2 * len);
    for (int i = 0; i < len; i++) s[i] = p[i] * 32768;
    write_file(name, s, len * 2);
    free(s);
}

#define FFT_SIZE 2048       // FFT 点数（建议为2的幂次）
#define SAMPLE_RATE 44100   // 音频采样率
int main() {
        // 1. 读取 PCM 数据（假设为16位有符号整型）
    float sample[FFT_SIZE];    
    // generate_triangular_wave(sample, FFT_SIZE, 16000);
    // generate_sine_wave(sample, FFT_SIZE, 16000);
    // write_float("input_data.pcm", sample, FFT_SIZE);
    FILE *fd = fopen("pcm_input.pcm", "r");
    for (int i = 0; i <FFT_SIZE; i++) {
        int16_t s;
        fread(&s, 1, 2, fd);
        sample[i] = (float)s / 32768;
    }
    write_file("input_data2.pcm", sample, 4 * FFT_SIZE);
    fclose(fd);
    
        // 2. 初始化 FFT 配置（实数FFT）
        kiss_fftr_cfg cfg = kiss_fftr_alloc(FFT_SIZE, 0, NULL, NULL);
    
        // 3. 输入和输出缓冲区
        kiss_fft_scalar *fft_in = (kiss_fft_scalar*)malloc(FFT_SIZE * sizeof(kiss_fft_scalar));
        kiss_fft_cpx *fft_out = (kiss_fft_cpx*)malloc((FFT_SIZE/2 + 1) * sizeof(kiss_fft_cpx));
    
        // 4. 分帧处理音频（示例：取第一帧）
        int frame_start = 0;
        for (int i = 0; i < FFT_SIZE; i++) {
            // 将 PCM 数据转换为浮点并归一化到 [-1, 1]
            fft_in[i] = (float)sample[frame_start + i];
    
            // 可选：应用汉宁窗减少频谱泄漏
            // float hanning = 0.5f * (1 - cos(2 * M_PI * i / (FFT_SIZE - 1)));
            // fft_in[i] *= hanning;
        }
    
        // 5. 执行 FFT
        kiss_fftr(cfg, fft_in, fft_out);
    
        // 6. 计算幅值频谱
        float spectrum[FFT_SIZE/2 + 1];
        for (int i = 0; i < FFT_SIZE/2 + 1; i++) {
            spectrum[i] = sqrtf(fft_out[i].r * fft_out[i].r + fft_out[i].i * fft_out[i].i);
            if (spectrum[i] > 1.02)
            printf("Spectrum  %dHz :%f\n",  i * SAMPLE_RATE / FFT_SIZE, spectrum[i]);
        }
    
        // 7. 清理资源
        free(cfg);
        free(fft_in);
        free(fft_out);
    
}